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ABSTRACT 



We have performed fully self-consistent A^-body simulations of star clusters near the Galactic 
center (GC). Such simulations have not been performed because it is difficult to perform fast 
and accurate simulations of such systems using conventional methods. We used the Bridge code, 
which integrates the parent galaxy using the tree algorithm and the star cluster using the fourth- 
order Hermite scheme with individual timestep. The interaction between the parent galaxy and 
the star cluster is calculate with the tree algorithm. Therefore, the Bridge code can handle both 
the orbital and internal evolutions of star clusters correctly at the same time. We investigated the 
evolution of star clusters using the Bridge code and compared the results with previous studies. 
We found that 1) the inspiral timescale of the star clusters is shorter than that obtained with 
"traditional" simulations, in which the orbital evolution of star clusters is calculated analytically 
using the dynamical friction formula and 2) the core collapse of the star cluster increases the 
core density and help the cluster survive. The initial conditions of star clusters is not so severe 
as previously suggested. 



Subject headings: galaxy: star clusters 
ical — stellar dynamics 

Introduction 



A few dozens of very young and massive stars 
have be en found in the central parsec of the 



Galaxy: center, kinematics and dynamics — methods: numer- 



the stars on the clockwis e rotating disk are circu- 
lar ( Paumard et al. 2006) or eccentric have lower- 
limit eccentricities of 0.0 — 0.8 ( Lu et a l. 2006), 



Gala xy (jKrabbe et al.lll995t iPaumard et al.ll2001 
20061 ). These stars are a few mi llion years old 



while those on the countercl ockwise disk have hig h 
eccentricities at around 0.8 (Pau mard et al . 2006). 
In the central parsec, in situ formation of these 
stars seems difficult because of the strong tidal 
field of the central BH. To overcome this difficulty, 
two possibilities have been suggested: (1) in situ 
star formation in a massive accretion disk, or (2) 
inspiraling young star clusters. 

The accretion disk scena rio was proposed by 
Levin fc Beloborodovi (|20031 ). A dense gaseous 



(IPaumard et~all l200ll iGhed |2003[) and lie on a 



disk ( Lu et al.ll2006l ) or two disks (|Paumard et al 
2006). The disks rotate around the central black 
hole (BH) and are at large angles with each other. 
One disk rotates cl ockwise in projection the other 
counterclockwise ( Paumard et al.l [2006). These 



disks are coeval to within 1 Myr. The orbit of 
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disk is formed from a molecular cloud which some- 
how fell to the neighborhood of the central BH. If 
the disk is sufficiently massive it can become gravi- 
tationally unstable, resulting in fragmentation and 
formation of stars. However, this scenario is prob- 
lematic. Observations have shown that two disks 
are at large angles with respect to each other and 
these stars on the disks formed almost at the same 
time. In this scenario, two disks must have ex- 
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isted simultaneously within 1 Myr. Moreover, it 
is difficult to m ake stars with eccentric orbit from 
accretion disks (Naya kshin et al 1 l2007h . 

The star cluster inspiral scenario was proposed 
bv lGerhardl (f^OOlh . A star cluster was formed at 
a distance of tens of parsecs from the GC and 
spiraled into the GC due to the dynamical fric- 
tion. This scenario is supported by the obser- 
vational fact that two young dense star clusters, 
the Arches and Quintuplet clusters are ob served 
at the distances of ~ 30 pc from the GC (|Figer 



2004) . Although this scenario can explain two stel- 
lar disks without difficulty, numerical simulations 
have shown that it would take too long time for 
the star cluster to inspiral to the central parsec 
unless it was very massive or its initial position 
is very near from the GC (Portegies Zwart et al. 
2003, hereafter PZ03; Gurkan & Rasio 2005). 

In these works (PZ03; Gurkan & Rasio 2005), 
the orbit of the cluster within its parent galaxy 
was c alculated using the d ynamical friction for- 
mula ( Chandrasekhar 19431 ). With this approach, 
the inspir al timescale might have been overesti- 
mated. In iFuiii et al.1 ( 2006f ). we performed fully 
self-consistent A-body simulations of a satellite 
galaxy within its parent galaxy and found that the 
orbital decay of the satellite is much faster than 
those calculated analytically from the dynamical 
friction formula. This difference was caused by 
particles escaped from the satellite. One mecha- 
nism is that the direct gravitational forces from 
escaped particles worked as effective drag force 
to the satellite. The second mechanism is that 
escaped particles remain close to the body of 
the satellite and enhance the dynamical friction. 
These effect should also work in the case of star 
clusters. Therefore, a fully self-consistent A-body 
simulation is necessary to obtain correct results 
for the orbital evolution of star clusters. 

Kim & Morris (2003; hereafter KM03) per- 
formed self-consistent A-body simulations of star 
clusters near the GC. Their results showed that if 
the initial central density of a star cluster was ini- 
tially very high (~ 10 8 M Q pc~ 3 ), the cluster can 
deliver stars to the central parsec of the Galaxy. 
In these simulations, however, the internal evolu- 
tion of the star clusters was neglected. The stars 
in their model of star clusters have an equal mass 
and a large softening length of 0.025pc. Such star 
clusters experience neither mass segregation nor 



core collapse. However, if the core collapse oc- 
curs, the core density of a star cluster increases. 
The initial high density of the core will not be 
necessary. Thus, it is also important to solve the 
internal evolution correctly. 

Such a fully self-consistent A-body simulation 
has been impossible with conventional numerical 
methods. While star clusters need a very accu- 
rate scheme such as the combination of fourth- 
order Hermite scheme and direct force calcula- 
tion, galaxies contain too many particles to use 
the direct force calculation. To solve this prob- 
lem, we have developed a new tree-direct hybri d 
scheme, the "Bridge" scheme (jFujii et al.ll2007t l. 
The Bridge scheme enables us to perform fully 
self-consistent A-body simulations of star clusters 
within their parent galaxies in a realistic time (less 
than 2 days with a single GRAPE-6 board). 

We performed fully self-consistent A-body sim- 
ulations of evolution of star clusters within their 
parent galaxies using the Bridge scheme. We also 
performed the "traditional" A-body simulations, 
in which the orbital decay of the star cluster is cal- 
culated using the dynamical friction formula for 
comparison. We found that the inspiral timescale 
of the star cluster is shorter than that obtained in 
previous studies. In addition, if the initial orbit of 
the star cluster is eccentric, the timescale of inspi- 
ral is much shorter than that for a cluster in the 
circular orbit of the same apocenter distance. 

The eccentricities of the stars escaped from the 
star cluster distribute around the eccentricity of 
the star cluster. Thus, if a star cluster is initially in 
an eccentric orbit, it naturally explains the rather 
high eccentricities of the stars in the observed 
"disks". Also, because of the mass segregation 
effect, very massive stars (more than 10M Q ) re- 
mained in the cluster and were brought very close 
to GC. Thus, we conclude that the timescale prob- 
lem with the cluster inspiral scenario was partly 
because of the wrong treatment of the dynamical 
friction and partly because the limited assump- 
tion of the circular orbit, and it is not difficult to 
make the central stellar disks from inspiraled star 
clusters. 

We describe the simulation method and initial 
conditions in section 2. In section 3 we show the 
results of simulations. Section 4 is for summary 
and discussions. 
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2. Numerical Simulation 
2.1. Models 

We adopted a King model with Wq = 3 for the 
model of a star cluster. Its core radius, r c , half- 
mass radius, r^, and tidal radius, r t , are 0.087 pc, 
0.13 pc, and 0.47 pc, respectively. It consists of 
65536 stars and we a ssigned each sta r a mass ran- 
domly drawn from a ISalpeterl (Il955h initial mass 
function between 0.3 and 100 M©, irrespective of 
position. The total mass of the star cluster, Mgc, 
is 7.9 x 1 4 MfTy Th is model imitates the Arches 



cluster (jNaeata et al.l 119951 ). whose mass and ve- 
locity dispersion are 7 x 1 4 M ( t i within 0.23 pc and 
22 km s" 1 , respectively (|Figer et al.l l2002h . The 



Arches cluster is located at ~ 30 pc from the GC. 
We, however, placed our cluster much closer to the 
GC, to compare our result with those of PZ03. In 
table [TJ we summarize the model parameters. 

We used two galaxy models as a model of the 
central region of the Galaxy; one includes the cen- 
tral SMBH (galaxy 1) and the other does not 
(galaxy 2). For the galaxy model 1, we adopted a 
King model with non-dimensional central poten- 
tial of Wo = 10. We scaled the density and veloc- 
ity dispersion of our model at 5 pc from the GC 
to the those of the Galaxy. The density at 5 pc is 
6.8 x 10 3 M©pc -3 and the one-dimensional velocity 
dispersion is 64 km s _1 in our scales, while the ob- 
served density at 5 p c is p(5pc) = 6.9 x 10 3 M©pc~ 3 
(jGenzel et al.ll2003l) and the o bserved velocity dis - 



persion is 54 km s 1 at ~ 4 pc (jGenzel et al . 2000). 
Figure [TJ shows the enclo sed masses of our m odel 
galaxy and the result of Genzel et al.l (2003). In 
this model, the total mass of the galaxy, Mq, is 
8 x 10 7 M Q , and the core radius, r Ci G, and the 
half-mass radius, r^Gi are 0.66 pc and 21 pc, re- 
spectively. We used 2 x 10 6 particles to model 
the Galaxy. The mass of a particle of the galaxy 
is 40 M Q . Since our galaxy model has the finite 
core size of 0.66pc, the orbital evolution of the star 
cluster should be reasonably accurate as far as its 
distance from the GC is more than lpc. 

The galaxy model 2 is a more realistic model of 
than model 1. It includes a central super-massive 
black hole (SMBH). This model is based on King 
model. We put a SMBH at the center of the galaxy 
with King model Wq = 10 and integrated it for 
around two crossing times. The central region of 
the galaxy evolved and its density profile became 



cuspy. This model roughly represents the Galactic 
center between 0.1 and 5 pc. The enclosed mass is 
shown in figure [TJ The total mass and the particle 
mass of the galaxy are 2.9 x 10 7 M Q and 14 M@, 
respectively. 

The initial position and velocity of the star clus- 
ter are shown in table [2 We calculated two orbits; 
one is circular (model C) and the other is eccentric 
(model E) for the galaxy model 1 , and two eccen- 
tric orbits (model Bl, B2) for the galaxy model 
2. In Model Bl, the star cluster has almost the 
same eccentricity as the previous simulation with- 
out SMBH. In Model B2, the eccentricity of the 
star cluster is lower than Model Bl. The orbital 
elements for the eccentric case were chosen so that 
the star cluster would survive at least for several 
orbits. At first, we tried the eccentric orbit with 
the same orbital energy as that of the circular case. 
However, in this case, the disruption of the cluster 
was too fast, unless the orbit is close to circular. 
Therefore we made the eccentric orbit significantly 
wider. 



The Galactic center 
Galaxy 1 
Galaxy 2 



0.1 



10 



R (pc) 



Fig. 1. — Enclosed mass of our model galaxy and 
the GC (jGenzel et al.j|2003f) . 
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Table 1: Models for the galaxy and the star cluster 

King W N M(M Q ) M bh (M ) r c (pc) r h (pc) r t (pc) 

Star cluster 3 65536 7.9 x 10 4 - 0~087 0T3 OAT 

Galaxy 1 10 2 x 10 6 8.0 x 10 7 - 0.66 21 120 

Galaxy 2 10 2 x 10 6 2.9 x 10 7 3.6 x 10 6 - 9.6 72 



Table 2: Initial Conditions for the cluster orbit 

Simulation Galaxy model Orbit Initial position (pc) Initial velocity (km/s) 
Model C I Circular 2 130 

Model E 1 Eccentric 5 72 

Model Bl 2 Eccentric 5 57 

Model B2 2 Eccentric 5 67 
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2.2. Fully Self-consistent iV-body Simula- 
tion 

We performed fully self-consi stent TV-body sim - 
ulations using the Bridge code (jFuiii et al.ll2007l h 
It is a direct-tree hybrid code. Only the inter- 
nal motion of the star cluster is calculated by the 
direct scheme with high accuracy, and all other 
interactions are calculated by the tree algorithm. 
The splitting between the direct part and tree part 
is through the splitting of the Hamiltonian in the 
way similar to MVS dWisdom fc Holmanl Il99lt 
iKinoshita. Yoshida. fc Nakail Il99lh or RESPA 
(jTuckerman et al.l Il990h . Thus, the low-accuracy 
calculation of galaxy particles and its interaction 
with cluster particles are integrated with time- 
symmetric leapfrog algorithm, resulting in small 
long-term error. Hence, we can treat a large-A" 
system with embedded small-scale systems fully 
self-consistently and accurately. 

The numerical parameters used for the time 
integration are summarized in table [3j For the 
tree part, the Bridge code needs the same pa- 
rameters with the Barnes-Hut tre ecode modified 
for the use with GR APE hardware (Bar nes fc Hull 
Il986t iMakinol 12004) . We used the opening angle 
9 = 0.75 with the center-of-mass (dipole-accurate) 
approximation. The maximum gr oup size for a 
GRAPE calculation (jMakinol Il99ll ) is 8192. The 
stepsize of leapfrog integrator is At — 1/512 (in 
Heggie unit) = 2.9 x 10" 4 (Myr) for galaxy 1 
and At = 1/1024 (in Heggie unit) = 1.2 x 10" 4 
(Myr) for galaxy 2. The potential is softened us- 
ing Plummer softening. The softening length be- 
tween galaxy particles that between cluster parti- 
cles and galaxy particles are the same. Both are 
6q = 3.9 x 10~ 2 pc. For galaxy 2, we adopted the 
softening length between between the SMBH and 
galaxy particles, £g-bh, is 0.12 pc and between 
the SMBH and star cluster particles, esc-BH, is 
0.012 pc. 

For the direct part, we used the fourth-order 
Hermite integrator with block timestep, and the 
timestep criterion is of th e standard Aarseth type 
(|Makino fc Aarsetbl[l992T ) with 77 = 0.01. We also 
used the Plummer softening for the gravitational 
force between cluster particles, and we did not 
model physical collisions or binary formation in 
the calculation reported in this paper. The soft- 
ening length between star cluster particles, esc, is 



1.0 x 10- 5 



pc. 



We stopped the simulations at T = 0.75 — 0.8 
(Myr). Since we used the softening and did not 
model the physical collision and merging of the 
stars, the structure of the star cluster after the 
core collapse might not expressed correctly. We 
ignored the stellar evolution because we treat only 
very short time (< 1 Myr) in our simulations. 

We used GRAPE-6 (|Makino et all l2003h for 
force calculation. The total energy was conserved 
better than 5 x 10~ 5 for the circular orbit 8 x 10~ 5 
for the eccentric orbit throughout the simulations. 

Table 3: Parameters for A-body Simulation 



Parameters 


Value 


ec 


3.9 x 10-* (pc) 


esc 


1.0 x 10~ 5 (pc) 
1.2 x 10" 1 (pc) 
1.2 x 10~ 2 (pc) 


CG-BH 


CSC-BH 


At (for galaxy 1) 
At (for galaxy 2) 


2.9 x 10" 4 (Myr) 
1.2 x 10" 4 (Myr) 


e 


0.75 




8192 



2.3. 



A-body Simulation with Artificial 
Dynamical Friction 



To compare our result with those of previous 
works, for model C and E, we also performed N- 
body simulations in which the Galaxy is modeled 
as a fixed potential and the dynamical friction 
due to the Galaxy is calculated analytically. This 
treatment is the same as in PZ03, where the accel- 
eration due to the dynamical friction is calculated 
using Chandrasekhar's dynamical friction formula 
dChandrasekharll943HBinney and Tremainelll987t 
iMcMillan fc Portegies ZwartJl2003l) 



2 ^sc 
a d f = -47rlnAxG pcMsc^— ■ 



(1) 



Here Msc and v§c are the mass and the center- 
of-mass velocity of the star cluster, pq is the local 
density of the Galaxy, In A is the Coulomb loga- 
rithm, and 



X = erf(X) 



2X 



exp(-A 2 ), 



(2) 
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where X = vsc/V^ a G an d &G is the local one- 
dimensional velocity dispersion of the Galaxy, as- 
sumed to be isotropic and locally Maxwellian. For 
Msc, we adopted the bound mass, which is the 
total mass of the particles bound to the star clus- 
ter. We gave all bound stars the acceleration due 
to the dynamical friction from the above formula. 
Unbound stars were not affected by the dynamical 
friction. 

In equation |T|), the Coulomb logarithm, In A, 
is given by 



In A = In 



(3) 



mum impact parameters. It is often set as a con- 
stant given by 



(4) 



where Rsc is the distance of the star cluster from 
the GC, (rsc) is the characteristic radius of the 
star cluster (roughly the half-mass radius). PZ03 
adopted a constant valu e, x m A = 1- However, 
Hashimoto et al. I (|2003l ) found a constant A over- 
estimates dynamical friction at pericenter and pro- 
posed a variable A: 



In A = In 



( Rsc \ 



(5) 



where esc is the size of the star cluster. We per- 
formed A^-body simulations with artificial dynam- 
ical friction in two ways. One is the same as 
that used PZ03 (c onstant A), and the othe r is the 
way proposed by lHashimoto et al.~ | (|2003h (vari- 
able A). We adopted the virial radius of the star 
cluster for esc- 

3. Simulation Results 



ring-like structures. Finally they form a disk-like 
structure. 

We investigated the orbital and internal evolu- 
tions of the star cluster. Figure [3] shows the dis- 
tance of the star cluster from the GC obtained by 
our A'-body simulations. The solid curve shows 
the result of the full A^-body simulation. The 
dashed and dotted curves show the result of the 
"traditional" simulations in which the dynami- 
cal friction is calculated analytically from the for- 
mula with constant A and variable A, respectively. 
The orbital decay in the full A^-body simulation is 
faster by 30-40% than that in other simulations. 
On the other hand, the evolution of the bound 
mass of the star cluster shown in figure 0] is almost 
the same among the three. This result suggests 
that previous studies underestimated the inspiral 
timescale of star clusters. This effect is not very 
large, but not negligible. 

In figure 21 it seems that the mass loss in the 
very late stage (after more than 80% of the ini- 
tial total mass is lost) seems to be significantly 
slower for the case of full A^-body simulation than 
those for other two simulations, even though the 
cluster is much closer to GC. This difference is 
probably due to the finite core size of our galaxy 
model (r c = 0.66pc) and not the reality. In order 
to study the evolution after this stage, we need to 
use a more realistic model of the mass distribution 
of the central parsec of the galaxy, which includes 
the central massive black hole. In addition, stellar 
collision and merger within the star cluster must 
be modeled. We are currently working on such an 
extension of the Bridge code. 

Figures [5] and [5] show the core radius and the 
core density of the star cluster. The core collapse 
occurred at around T ~ 0.53 Myr. The core col- 
lapse times are also the same among the three sim- 
ulations. 



3.1. Circular Orbits 

The top panel of figure [2] shows the snapshots of 
the star clusters projected onto the x-y plane. The 
top and bottom panels show model C (circular or- 
bit) and E (eccentric orbit), respectively. In model 
C, the star cluster is initially located at 2 pc from 
the GC. Due to the tidal field of the Galaxy, the 
star cluster becomes elongated. Particles stripped 
from the star cluster form tidal arms and make 
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Fig. 2. — Snapshots of the star clusters projected onto x — y plane. The upper six panels are for the run 
with the circular initial orbit (model C), and the lower six panels are for the run with the eccentric orbit 
(model E). Times are 0.122, 0.239, 0.356, 0.473, 0.5917, and 0.708 Myrs. 



3.2. Eccentric Orbits 



The bottom panel of figure [2] shows the snap- 
shots of model E projected onto the x-y plane. The 
initial distance of the star cluster from the GC is 
5 pc. The star cluster is elongated and particles 
are stripped due to the tidal force of the Galaxy. 
The stripped particles form complex tidal tails. 

Figure [7] shows the orbital evolution of the star 
cluster. The orbital decay of the full iV-body sim- 
ulation is faster than the traditional simulations, 
as was the case in the runs from the circular orbit. 
In this case also, the evolution of the bound mass 
is the same among the three simulations (see fig- 
ure [5]). Moreover, this result shows that variable 
A works better than constant A. 

Figures [9] and [10] show the core radius and the 
core density of the star cluster. These show that 
the core collapse occurred at around T ~ 0.55 
Myr. The core collapse time is almost the same as 
that for the case of the circular orbit, even though 
the bound mass at the collapse time is different 
by almost a factor of two. The half-mass relax- 
ation time at the time of the collapse was 0.11 
Myr and 0.34 Myr, in model C and E, respec- 
tively. Thus, if clusters are rapidly losing its mass 
due to the tidal field, the apparent age measured 
by the present relaxation time can show large vari- 
ations, even if they started from the same initial 
condition and collapsed at the same time. Fur- 
thermore, the core collapse times are the same as 
that in the case without massloss due to the tidal 
field. The core collapse time, t cc , is estimated as 
t cc ~ 0.20t r h, where t r h is half-mass re laxation 
time (|Portegies Zwart k, McMillaiJ liooih . From 



this equation, we obtain the core collapse time of 
our model star cluster as 0.51 Myr. Our results 
agreed with it. 

Figure [11] - [14] shows the results of model Bl 
(dashed curves) and B2 (dotted curves), where the 
Galaxy model has a central BH. Solid curves show 
the result of model E (no BH; the same as the 
solid curve in figure 171- [TO]). The orbital evolution 
of the star cluster in model Bl and B2 is essentially 
similar to that in model E. However, the orbital 
decay in model Bl is somewhat slower than that 
in model E. As is clear in figure [T2l in model Bl, 
the mass loss at the pericenter passage is much 
larger than in the case of model E. Therefore, the 
core collapse did not occur and the core density 
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Variable A 
Constant A 
Full Af-body 



0.2 0.4 0.6 0.8 

T (Myr) 

Fig. 3. — The distance of the star cluster from the 
GC plotted as a function of time for model C. Solid 
curve shows the result of the full iV-body simula- 
tion. Dashed and dotted curves show the results of 
the "traditional" simulations with variable A and 
constant A, respectively. 




0.2 0.4 0.6 0.8 

T (Myr) 

Fig. 4. — The bound mass of the star cluster plot- 
ted as a function of time for model C. Curves have 
the same meanings as in figure [3] 
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Fig. 5. — The core radius of the star cluster, r c , 
plotted as a function of time for model C. Curves 
have the same meanings as in figure [3J 
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Fig. 7. — The distance of the star cluster from 
the GC plotted as a function of time for model 
E. Solid curve shows the result of the full TV-body 
simulation. Curves have the same meanings as in 
figure [3J 
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T (Myr) 



Fig. 6. — The core density of the star cluster, p c , 
plotted as a function of time for model C. Curves 
have the same meanings as in figure El 




0.2 0.4 0.6 0.8 



T (Myr) 

Fig. 8. — The bound mass of the star cluster plot- 
ted as a function of time for model E. Curves have 
the same meanings as in figure [7] 
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T (Myr) 

Fig. 9. — The core radius of the star cluster, r c , 
plotted as a function of time for model E. Curves 
have the same meanings as in figure El 




0.2 0.4 0.6 

T (Myr) 



0.8 



Fig. 10. — The core density of the star cluster, p c , 
plotted as a function of time for model E. Curves 
have the same meanings as in figure [71 



of the cluster did not increase. This is because 
the star cluster suffers the strong tidal force from 
the SMBH at the pericenter. The small mass of 
the star cluster slows the orbital evolution of the 
star cluster. Furthermore, the strong tidal field 
disrupts the star cluster much faster than the case 
without SMBH. 

In model B2, the pericenter of the star cluster 
is farther than that in model Bl. The evolution 
of the bound mass is the same as the case with- 
out SMBH. However, the orbital evolution is much 
slower than that of model Bl, because the pericen- 
ter of the star cluster is farther. 

The cluster model and its orbital evolutions in 
model E and Bl are almost the same as simula- 
tion 8 of KM03. Their model of the Galaxy has a 
power-law density profile with the central SMBH, 
while our galaxy model 1 for model E has no 
SMBH. However, the enclosed mass of the galaxy 
in KM03 is very similar to ours between around 1 
and lOpc. In their simulation, the star cluster was 
totally disrupted before 1 Myr, but in our simu- 
lation, the star cluster has 30% of its initial mass 
at the end of our simulation (0.75 Myr). The dif- 
ference is caused by the internal evolution of the 
star cluster. In our simulation, mass segregation 
and core collapse of the cluster made the central 
density much higher, which prevent the complete 
disruption of the cluster. On the other hand, in 
KM03, such an evolution was prevented by the 
numerical method they used. 

On the other hand, the cluster in model Bl was 
disrupted before the core collapse occurs by the 
strong tidal field of the central SMBH. At around 
1 pc, the enclosed mass of galaxy model 2 which 
includes a SMBH is twice as large as that of K04's 
model. Such a difference of the enclosed masses 
is caused by the mass of the SMBH . We adopted 
3.6 x 1O 6 M (|Eisenhauer et al.ll2005h . while KM03 
did 2.5 x 1O 6 M . 

Thus, the presence of a SMBH has consider- 
able effect on the evolution of star clusters. The 
main effect is that the tidal field near the GC be- 
comes much stronger, resulting in faster mass loss 
from star clusters and preventing its core collapse. 
In the models we tested, the core collapse time is 
about the same as the time of the complete dis- 
ruption in model Bl. If the core collapse time 
is somewhat faster, it is possible that the cluster 
survives and approaches to the GC. In particular, 
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if an IMBH is formed within the star cluster, it 
would significantly help the survival of the clus- 
ter. We will study this aspect in more details in 
the forthcoming papers. 

3.3. Eccentricities and Inclinations of the 
Escaped Stars 

For some of the stars in the central parsec, the 
projected positions and the proper motions have 
been measured, and their orbital elements have 
been estimated ( Paumard et al. 2006t Lu et al 
20061 ). Young and bright stars apparently be- 



long to one of the two "disks" (clockwise and 
counter-clockwise rotating disks), though the exis- 
tence of th e cou nter-clockwise disk is controversial 
20061 ). The eccentricities of the stars 



(Lu et al 



on the counter-clockwise rot ating disk are high, 
" 120061) . For the stars on 



0.8 (Paumard et al 



the clockwise rotating disk, Paumard et all (120061) 



concluded their orbits are circular, while Lu et al 



( 20061 ) concluded the lower limits of their eccen- 
tricities distribute between 0.0 and 0.8. These 
high eccentricities are difficult to explain with the 
in-situ formation scenario, and also have been 
thought to be difficult to explain with cluster in- 
spiral scenario, since in both cases the stars would 
have close-to-circular orbits. 

We investigated the eccentricities, e, and incli- 
nations, i, of stars escaped from the star cluster 
(i.e. unbound stars) for model C and E. Figures 
[TBI and [16] show the eccentricities and inclinations 
of escaped stars and their evolutions for model C 
and E, respectively. In these figures, the x-axis 
is the semi-major axis, a, of the star. The ec- 
centricity is defined as e = (r a — r p )/(r a + r p ), 
where r a and rp are the apocenter and percen- 
ter distances, respectively. We obtained the peri- 
and apocenter distances by integrating the orbits 
of the stars in the model potential. The inclina- 
tion is defined as i — cos^ 1 (h z /h), where h and h z 
are the angular momentum and its z component, 
respectively. The top panels show the eccentric- 
ities and inclinations of the stars escaped before 
T = 0.15 Myr. The central gap corresponds to 
the semi-major axis of the star cluster, and the 
left and right wings are the stars on the leading 
and trailing arms, respectively. The distributions 
of the escaped stars on the a — e and a — i plains 
expand because of the time-varying gravitational 
force from the star cluster (see the middle panels) . 
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Fig. 11. — The distance of the star cluster from 
the GC plotted as a function of time. Solid curve 
shows the result of the galaxy model without the 
central SMBH (Model E; the same as the solid 
curve in figure [7J). Dashed and dotted curves show 
the result of the galaxy model with the central 
SMBH (Model Bl and B2). 
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Fig. 12. — The bound mass of the star cluster 
plotted as a function of time for the galaxy model 
with the SMBH. Curves have the same meanings 
as in figure [TTJ 
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Fig. 13. — The core radius of the star cluster, r c , 
plotted as a function of time for the galaxy model 
with the SMBH. Curves have the same meanings 
as in figure [HJ 
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Fig. 14. — The core density of the star cluster, p c , 
plotted as a function of time for the galaxy model 
with the SMBH. Curves have the same meanings 
as in figure [Til 



The bottom panels show the distribution of all es- 
caped stars at the last snap shots, T = 0.75 Myr. 
Even in the case of the circular orbit (model C), 
the eccentricities of the escapers in innermost or- 
bits can reach rather large values. However, in the 
case of the eccentric orbit (model E), the eccen- 
tricities of the escapers are even higher. They dis- 
tribute between 0.4 and 0.8, being consistent with 
the observed values of stars in the two "disks". 
Thus, if the cluster is initially in a highly eccen- 
tric orbit, the high eccentricities of observed stars 
are naturally explained. 
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Fig. 15. — The eccentricity (left) and inclination (right) of the stars escaped from the star cluster as the 
function of their the semi-major axis a for model C (circular orbit). The top panels show the unbound stars 
at T = 0.15 Myr and the middle panels show the position of the same stars shown in top panels, but at 
T = 0.75 Myr. The bottom panels show the all unbound stars at T = 0.75 Myr. 
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Fig. 16. — Same as figure [T5l but for model E (eccentric orbit). 
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3.4. The Evolution of the Mass Function 



rewritten as 



Figure [T7l shows the evolution of the mass func- 
tion (MF) of the stars bound to the star cluster for 
model C. We can see that the slope of the mass 
function for to > 1OM becomes flatter as the 
system evolves, while that for to < 1OM shows 
rather little change. To see this effect more clearly, 
in figure [TH] we plot the fraction of the mass re- 
tained in the cluster as the function of the stellar 
mass for model C. Stars with mass to > 30M© 
are almost perfectly retained in the cluster, and 
the retention rate quickly drops in the range of 
10M© < m < 30M S . For m < 10M Q , retention 
rate becomes smaller for smaller mass, but the de- 
pendence becomes much weaker. 

Figure [19] shows the evolution of the enclosed 
mass in r, the distance from the center of the star 
cluster for model C. The enclosed masses are cal- 
culated for five mass ranges. Initially, all mass 
ranges have same profile. At T — 0.15 Myr, the 
most massive stars have sank to the center and the 
second massive stars also shows some central con- 
densation. At T = 0.45 Myr, a large fraction of 
stars more massive than 10 Mq have sank to the 
center (within radius 0.02 pc). However, the dis- 
tribution of stars with mass less than lOM© have 
not changed significantly. Thus, the stars heavier 
than 10 Mq remained in the star cluster. 

We can estimate the critical mass of the star at 
which this change in the behavior occurs, by cal- 
culating the mass of the star at which the dynam- 
ical friction just balances the two-body heating. 
First, we consider the energy change of the stars in 
the system consisting of two components with the 
masses mi and m 2 - Assuming that the velocity 
dispersions of both components are Maxwellian, 
the mean energy change of a star of mass mi is 
expressed as 



- {miEl ) 



4v37rG mi 7712^2 log A 



-(m lEl ) 



An2m\mi(mn — mi), 



(7) 



where A is a constant. Now we consider the case 
of continuous mass distribution, in which the mass 
distribution of mi is given by 



dn 
dm 



(8) 



where C is a constant. The number density of 
mass m,2, 112, is expressed as 



ri2 = Cm 2 



(9) 



By substituting equation §§§ to ([7]) and integrating 
it over mass m2, we can obtain the energy change 
of stars with mass mi as 



df 



(mi£7i) 



A' mi 



cc mi 



1 / 




J m m 




— OL- 


h3 






-a + 


3 



m^ a+1 (m 2 - mi)^ 



-mi 



(U) 



where m max and m m i n are the maximum and min- 
imum mass of the MF and A' is a constant. If the 
right side of equation (fTTj) is negative, the star 
with mass mi loses energy and sink to the center 
of the star cluster. The minimum mass with the 
negative energy change, m S i n k, is expressed as 



T^sink 



/(a)m„ 



m~ a+2 

/,A max 



TO, 



-a+2 



1 - 2-- tt + 3 

1 -x- a + 2 ' 



(12) 
(13) 



(to 2 (E 2 ) - mi (E 



where we defined x = m m i n /m max and f(a) = 
(—a + 2)/(— a + 3). Using the values used for our 
model, m max = 100M Q , m m i n = O.3M and a — 
2.35, we obtain TO s j n k = 7.9M . This value agrees 
;cll with the mass where the power-low index of 
e MF breaks in figure [TH] If x -C 1, we obtain 



where (Ei) and (E2) are mean specific kinetic en- 
ergy of each components, n 2 is the number den- 
sity of the stars o f mass to 2 , and In A is the 
Coulomb logarithm (jHeggie fc H ut 2003). At least 



in the initial model, the velocity dispersion of the 
stars is independent of the mass. So we can set 
(Ei) = (E 2 ) = 1 /2a 2 . In this case, equation is 



1 



-a+3 



r-a+2 



Therefore, 



m-sink 




1 (a > 3), 

a - 2 (2<a< 3), (14) 
(a < 2). 



f(a)m 

-f(a)m 

f(ot)m 



a-2 3-q 
min max 



(a > 3), 

(2 < a < 3C3,5) 

(a < 2). 
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When 2 < a < 3, the value of m s i n k varies from 

^min to /71 max . 

The observed MF of the stars in the central par- 
sec is much flatter than Salpeter ( Paumard et al.l 
2006) . The cluster inspiral model rather naturally 
explain this flat MF, since only the most massive 
stars remain bound to the cluster and are carried 
to the central region of the galaxy. 

4. Summary and Discussion 
4.1. Summary 

We performed fully self-consistent iV-body sim- 
ulations of a star cluster within its parent galaxy 
and compared the orbital and internal evolutions 
of the star cluster with those obtained by "tradi- 
tional" simulations, in which the orbital evolution 
of the star cluster is calculated from the dynam- 
ical friction formula. We confirmed that the in- 
spiral timescale of the star cluster is shorter than 
that obtained from the "traditional" simulations. 
Furthermore, our results showed that the core col- 
lapse make the core density of the cluster increase 
and helps the cluster survive. 

We performed simulations of circular and eccen- 
tric orbits of the star cluster. In previous studies, 
most of the simulations were from circular orbits 
(PZ03; Giirkan & Rasio 2005). We found that, 
however, eccentric orbits are more favorable to ex- 
plain the distribution of stars around the GC for 
following reasons. 

First, eccentric orbits are natural, if the for- 
mation of star clusters was triggered by collisions 
between gas clouds. Second, star clusters with ec- 
centric orbits can approach to the GC much faster 
than those wit h circu lar orbits (KM03). Third, 
Paumard et al. I (|2006h showed that many stars in 
the counterclockwise rotating disk have high ec- 
centricities (e ~ 0.8), while the distribution of the 
eccentricities of the stars in the clockwise disk is 
very broad. Since the eccentricities of the escaped 
stars distribute around the eccentricity of the star 
cluster, the star cluster model with eccentric or- 
bits can naturally explain the existence of high- 
eccentricity stars. 

The power-low index of the MF of the bound 
stars to the star cluster breaks at around 7.9 Mq. 
Stars heavier than this mass sink to the center of 
the star cluster due to the mass segregation. Since 
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Fig. 17. — The mass function of the bound stars 
for model C. Times are T = 0.15, 0.3, 0.45, 0.6, 
and 0.75 Myrs. 
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Fig. 18. — The fraction of the mass function of the 
bound stars to the initial mass function for model 
C. Times are the same as figure [j~7l 
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Fig. 19. — The evolution of the enclosed mass of 
the bound stars for model C. 



the tidal stripping removes the stars outside of the 
star cluster, the massive stars selectively remain in 
the star cluster. As a result, the star cluster car- 
ries only massive stars to the GC. The star cluster 
scenario can reproduce the flat MF in the central 
parsec, without the need for nonstandard initial 
mass function. 

4.2. Realistic Model of the Galaxy 

First, we showed the orbits of star clusters de- 
cays faster in full iV-body simulations than in "tra- 
ditional" simulations, which treat the dynamical 
friction analytically using a King model Wo = 10 
as a model of the central region of the Galaxy. The 
model is sufficient for such comparisons, but not 
for more realistic comparison between our model 
and the stars in the GC because within ~1 pc 
the mass density of our model is much lower than 
that of the actual Galaxy. Next, we showed the 
case of more realistic galaxy model with a central 
BH. The orbital evolution was similar that in the 
model without a BH. For the comparison between 
our simulations and the actual stars, however, we 
need more simulations in various initial conditions. 
We will report more detail result of runs with the 
central SMBH in forthcoming papers. 

4.3. Formation of an Intermediate-mass 
Black Hole (IMBH) 

IRS 13E consists of seven stars within a pro- 
jected diameter of ~ 0.02 pc and is located at ~ 
0.14 pc in projection from the GC and these stars 



have very similar proper motions ( | Maillard et al 
2004t|Paumard et al.ll2006h . iMaillard et all \2004) 



suggested that IRS 13E is the remnant core of a 
star cluster that have fallen to the GC and dis- 
solved there and that the members of IRS 13E are 
bound by a central IMBH. From the analysis of 
proper motions, the minimu m mass of the IMBH 
was estim ated as 1300 (IMaillard et all [20041 ) 



1O 4 M (ISchodel et al.ll2005h 



In this paper we simulated the evolution of 
star clusters only before the core collapse because 
of the limitation of our present code. Our code 
currently cannot treat the post-collapse evolution 
since we use the softened potential. Collisions and 
mergers between stars would have occurred and 
an IMBH would have be formed, if our code can 
handle these events. The star cluster inspiral scc- 
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nario might reveal the origin of the IRS 13E. We 
are currently working to implement collisions and 
mergers. The result will be reported in the future 
papers. 
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